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INTRODUCTION 

Kamerlingh Onnes discovered some ninety years ago 
that the resistivity of certain metals would drop to zero 
when the temperature was lowered below a critical value 
Tc- He accordingly termed superconductors all materials 
which exhibited such state of disipationless electric flow. 
Superconductors display as a matter of fact a rich vari- 
ety of physical phenomena other than superflow, some of 
which are closely related to it, like the Meissner effect, 
some others having a different nature, as the electronic 
gap which appears in tunneling experiments. 

Superflow and Meissner's effects are a consequence of 
the London equation, a linear relation between the den- 
sity current and the vector potential, ] + KA = Q, which 
serves to define the current response tensor K. Such an 
equation, or its non-local version, holds in a supercon- 
ductor even in the presence of the usual scattering mech- 
anisms which degrade the current in a normal metal. 
A conventional exercise consists in solving London and 
Maxwell's equations in a specific geometry to find tiie. 
distribution of current densities and magnetic fields .Id'EI 
It is shown in this way that they are highly nonuniform 
across the sample: they are confined to a thin shell in 
contact to the surface and do not penetrate the bulk of 
the sample. 

Interestingly enough, superflow and Meissner's effects 
might be the most difficult properties to extract from 
the microscopic BCS theory of superconductivity, be- 
cause the derivation of the London equation is not an 
easy task. The standard procedure to derive this rela- 
tion uses linear response theory.l3 One complication arises 
from the fact that such expression is not gauge invariant, 
so that a change in the longitudinal part of A will change 
J/. The longitudinal response kernel K; is in turn related 
to the density response function through the continuity 
equation, so the problem of gauge invariance can also 
be stated as a problem of charge conservation. Within 



linear response theory, BCS theory must therefore be ex- 
tended to enforce the Ward Identity which relates the 
dressed vertex of theJongitudinal response function K/ to 
the BCS self-energyE It can also be proven that any ex- 
plicit computation of superflow in a superconductor must 
be fully seljE-jConsistent in order to ensure that charge is 
conserved. Qu 

The superflow state can be destroyed by increasing the 
current density, the magnetic field, or the temperature 
above certain critical values which depend among other 
things on the specific material used, the amount of dis- 
order of the sample, or the geometry of the device which 
is used. Furthermore, the relation between j* and A is 
not linear for large j*. In other words, London's equation 
does not hold close to the critical value of the super- 
fiow current density, jc. The situation in strong type II 
superconductors is even more complicated because the 
magnetic field originated by the supercurrent might give 
rise to vortices; these vortices movonin the presence of the 
supercurrent and dissipate energyfl 

Superflow in High-temperature superconductors might 
be anisotropic because (a) the order parameter has d- 
wave symmetry; (b) the electronic structure is also fairly 
anisotropic, witli van Hove singularities slightly below the 
Fermi energy;Ln3 (c) the two-dimensional Cu O2 planes, 
responsible for much of the physics of the cuprates, 
present different types of orthorhombic distortions. The 
purpose of this article is to study the physics of super- 
flow in clean, infinite, d-wave superconductors on an or- 
thorhombically distorted square lattice. A phenomeno- 
logical two-dimensional tight-binding model which cap- 
tures some of the experimental features of the supercon- 
ducting state of High-tcmpcrature superconductors will 
be used.EHld No attention will be paid to the effects of 
disorder, twin boundaries, or external or induced mag- 
netic fields, all of which can be important in this com- 
pounds and will be dealt with in the near future. Super- 
flow in two-dimensional d-wave superconductors has al- 
ready been studied using a Ginsburg-Landau approach, 
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which focused in particular on the induction of non-^L 
wave components upon switching on of the current .til 
Superflow in quasi-one-dimensional s-wave-suaerconduc- 
tors has also been stU|died| in equilibriumtjO and non- 
equilibrium situations .Mm 

Section II is devoted to present the model and the for- 
malism used to solve it in the presence of a finite super- 
flow. The mapping onto a linearized isotropic model is 
also shown. Section III discusses the results obtained for 
the phase diagram of the model, the gaps, their distor- 
tion due to orthorhombic effects and the local density of 
states, for zero current density. Section IV is devoted to 
study the doping and angular dependence of the gaps, 
free energy, current density, and critical current in the 
presence of a superflow. A summary closes the article. 
Even though the analytical calculations have been per- 
formed using the Matsubara formalism, all the numerical 
results provided in the article correspond to T=0. 



FORMALISM 

The partition function of a superconducting sample in 
the presence of an electromagnetic field can be written 
as 



Z[T,ii) = / P[A*,A,i*,Ao]e- 



-iS[A-,A,i*,ylo] 



(1) 



where the functional integral sums over all possible 
configurations of order parameter and electromagnetic 
fields, A(x,x') = |A(x,a;')|e'^'^(^)+'*(^'», Afx), Ao(a;), 
and where the imaginary-time action is 



dx dx' 



-Tr In G-\x,x') + 



'M2 



\A{x,x') 
V{x,x') 



(2) 

+ [enao{x) + ^^^^^f^^)six-x')\ +7^*(0). 



Here x denotes the four vector (r, x), V{x,x') = V{x — 
x') 5{t—t') is the pairing potential, e rii is the background 
charge density of ions, Ge{x,x') is the Green function of 
electrons in Nambu's space, 

(ndr+H(x)'^S{x-x') A{x,x')S{t-t') 
A*{x,x')5{t ^t') (hdr-n*ix)^ S{x-x') 

and in the Hamiltonian 

^ (p + e/cA(x))'^ , . . , , 

2 m 

the Zeeman term has been discarded. 

It can be seen that the action is invariant under the 
following gauge U(l) symmetry 
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Ao{x) ^ Aa{x) - —dr9{x), 
2 e 

A{x) ^ A{x) + —Ve{x), 
2 e 



(3) 



(x) i~> 4'{x) 



e{x)_ 
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up to a total time derivative. 

The phase of the order parameter field can there- 
fore be absorbed by the electromagnetic field, and drops 
out of the-.problem if one defines the gauge invariant 
quantitiesE^ 



1 fi 

*,(.t) = An{x) -—drHx), 

2 e 



Vs{x) 
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V0(a;) 
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Aix] 



(4) 



in terms of which the electric and magnetic fields are 
written as 



E{x) = -V*,(.t) drVsi^), 

e 



B{x) 



V X Vs{x). 



(5) 



As shown by van Otterlo and coworkers the whole 
problem can be rephrased using these new variables: 

Z{T,fi)= J P[A,,i/„vI/,]e-*'S'[A,A,*,]^ 

Afix,x') 



S' — J '^x' 
+ (eni'<ifs{x) - 



-Tr In G-^{x,x') + 



V{x,x') 



B^x)-E^{x) 
8^ 



S{x-x') +H*{0), 
Ai5{t 



(hdr + 1-1^ S{x - x') /\id[T-r 

AiS{t-t') (ndr-H*)6{x-x') 



^ 2- = ^ „ ^ " + Viattix) - e - fl, (6) 

2 m 

where A; = |A| is a real scalar field whose vacuum value 
will be different from zero below Tq and the phase (p is 
^ explicitly linked to the dynamics of the electromagnetic 
fields and the Higgs mechanism through the longitudinal 
, part of the superfluid velocity field, Vs,i- 

Stationarity of the action with respect to A; , Vs and 
leads to the saddle point equations which determine 
the classical configurations of electromagnetic and order 
parameter fields, and to the saddle point approximation 
for the thermodynamic potential 



V ■ Eci{x) = 47re(ni - ne{x)), 



(7a) 



- -drE,i{x) + V X B,i{x) = —Six), (7b) 
c c 
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Ai^ci{x,x') = -j-^V{x,x')Ge,i2{x,x'), (7c) 



^sp{T, fj,, Ai ci, Eel, Bel) = -wT^Sp[^i,ch Eci,Bci], (7d) 

p n 

where L'^ is the volume of the sample, and Ge,i2 is the 
ofl-diagonal element of Ge- The electronic and current 
densities are defined in terms of the diagonal elements of 
Go as 



ne{x) 



Ge,ii(a;,a;) - Ge,22(a;,a;) 



:Tr 



{pi + mvs{x)f:,) Ge{x,y) 



(8a) 
(8b) 



where fz is the third Pauli matrix. 

The saddle point equation for the order parameter can 
be further simplified by changing to center of mass and 
relative coordinates and performing the Fourier trans- 
form with respect to the relative coordinate. Then 

Ai,ci{k,R) = -^Yl V{k~k')GeMk',R)- (9) 



Expanding now Ai^ci and V in terms of some basis set 
of functions rja as can be, for instance, those which form 
the irreducible representations of a crystallographic point 
group. 



Ai^ci{k,R)=Y,^a{R) Va{k), 

a 

V{k ~ k', R)^J2 ^''(^) '?-(^) ^(^"') 

a 

one easily finds that 



AaiR) 



VajR) 



(10) 
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Superfiow is a superconducting state characterized by 
the stationary flow of a current j* in the absence of any 
longitudinal electric field (neglecting that built up at 
the boundary of the specimen.) This implies in view 
of Eqs. (^ and (|^) that (x) is identically zero, and 
that the superfiuid velocity Vs{x) does not depend on 
time. This state is therefore determined by two classical 
fields, Ai^ciix) and Vs^ciix), which must be determined 
self-consistently with the aid of Eqs. ( [7b| ) and (^). 

In an infinite sample, with no externally applied fields, 
Ai^ci and Vs^ci are uniform, so that Eq. ( [fh] ) drops out 
of the problem and the superfiuid velocity is just an in- 
put parameter, which can be viewed as a Doppler shift 
acquired by the kinetic energy part of the Hamiltonian 



r = 



(p + rnvs)^ _ {p + hq)'^ 
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The order parameter field must still be determined 
and, once self-consistency is achieved, the current density 
can be computed via Eq. (8b). The continuity equation, 
V • is automatically satisfied in this case. 

In a finite geometry with no constrictions, like a wire or 
a thin film, the full self-consistent problem must be han- 
dled, but there are some facts which simplify the compu- 
tation: (a) ^(x) must always point in the same direction, 
say the Z axis, because of the symmetry of the prob- 
lem, (b) The continuity equation then demands that J 
and Vs vary only in the XY plane. Now, the continuity 
equation is automatically satisfied if the order parameter 
is calculated self-consistently.D'EI Therefore, any explicit 
self-consistent calculation conserves charge and is gauge 
invariant by construction if j* is computed using Eq. ( pq ) 
with no further approximations like linear response the- 
ory. 

For an infinite superconductor, the current density can 
also be obtained in a more direct fashion by calculating 
the gradient of the thermodynamic potential. 



(13) 



It is worth pointing out that a finite Doppler shift in 
the dispersion relation of a normal metal leads to no ob- 
servable effects, because it can always be made disappear 
by performing a simple change of variables in the k-sums, 
k + q <—> k. Indeed, a finite current can only be induced 
in a metal by applying a finite electric field. 

The formalism developed so far is going to be ap- 
plied to study the physics of superfiow in a clean, in- 
finite d-wave superconductor where the electromagnetic 
fields are zero or at least negligible. The electronic 
fields will be described by a phenomcnological one-band 
Hubbard-like [toodel on an orthorhombically distorted 
square lattice .El'E^ 

The dispersion relation used comes from a fit to the 
spectral peaks found in AIIPES experiments of optimally 
doped Bi2Sr2CaCu2O8-5.0 The interactions among elec- 
trons V{x — x') are assumed to be different from zero only 
when X and x' are the same site {V(0) = Vq) or when they 
are nearest neighbors {V{x — x') = V{S) = Vi). Explic- 
itly, 



ct cr 

k.a K,' 



(14) 



where z is the coordination number, the sum in 6 extends 
up to five nearest neighbors, the hopping integrals t{\5\) 
have the values [-595,164,-52,56,51] meV, and the chemi- 
cal potential can be varied freely to change the mean oc- 
cupation number n of the band. It is therefore assumed 
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that the band is shifted rigidly upon doping, a fact which 
is not inconsistent with the available data from ARPES 
experiments.113 

Monte Carlo and RPA studies of the repulsive Hub- 
bard model have found that the effective potential among 
dressed quasiparticles assumes an RKKY or "Mexican 
hat" fornJl3 so that it has a large repulsive on-site con- 
tribution and then decays exponentially in an oscillating 
way. Within this_magnetic scenario for superconductivity 
in the cupratesjlj the length scale over which the poten- 
tial varies is set by the magnetic coherence length S^mag 
and is therefore uncommensurated with the lattice spac- 
ing a. The interaction term used in this article bears 
some resemblance with such an effective potential if the 
on-site interaction Vo if assumed to be large and negative. 
The resulting effective model can then be thought of as 
a low energy fixed point Hamiltonian which appears af- 
ter integrating out the high frequency degrees of freedom 
of a more complicated model describing the physics of 
the cuprates. The effective model can therefore be safely 
solved using Mean Field techniques. 

Generically, the four possible one-dimensional irre- 
ducible representations of the group of the square C4t, 
transform under its symmetry reoperations as the identity, 

— ^ xy and xy {x"^ — y^)." Their explicit form in a 
tight-binding formulation would be s or extended s (s*), 
dg.2_y2, dxy, and dxy x dx2_y2, states. 

A Mean Field decomposition of the on-site interaction 
naturally gives rise to a pure s or to a spin-density-wave 
(SDW) state, depending on whether Vb is attractive or 
repulsive. A nearest neighbor attraction produces an s* 
and a dx2-y2 state; a second nearest neighbors interac- 
tion V2 gives rise to a second s* state and to the dxy state, 
and so on. The Mean Field solution of this Hamiltonian 
therefore allows for a SDW state, and superconducting 
s-wave, extended s-wave, d-wave states, and mixed rep- 
resentations of them [s + s* ^d + i{s + s*), s + s* + d). 
Similar results have been obtained in Refs. ^-p^. 

The saddle point equations for the thermodynamic po- 
tential, the order parameter, the occupation of the band 
and the Helmholtz free energy per site are 



^lsp{T,^i) I A, 



M 



I A, 



-E 

k 



Vi 



^k+q 



k,q 



■ln(l 



M ^ 2Er - 

t k,q 



tanh 



N 



k L 



1- ^^±^-^tanh'^^^^'« 



^k+n + ^k 



2Ei: - 

k.q 



V 2 



(15a) 

(15b) 
, (15c) 



fsp[T,n) = — \-^J.{n)n, (15d) 



TABLE I. Range of variation of the interactions along with 
their respective order parameters and 77-functions 



Parameters 


Range (meV) 


Gaps 






Vo 


[-400,400] 


As 


1 




Vi 


[0,1500] 


As, 


cos(kx£,) + 


cos{kyi) 


Vi 


[0,1500] 


Ad 


cos(kx^) - 


cos{ky^) 



where the dispersion relation for the quasiparticles is 



Er.. 




(16) 



a— s,s* .d 



and the 77-functions, which form a truncated tight- 
binding basis of the crystal point group C^y, are given 
in table I, along with the allowed range of variation of Vb 
and Vi. It has been assumed in the numerical calcula- 
tions of this article that ^ is equal to the lattice spacing a, 
but it can generically be any other characteristic length 
scale, like the magnetic coherence length S,mag- 

Orthorhombic effects are present in most of the 
cuprates, giving rise to a mixing of the irreducible rep- 
resentations. For instance, in YBa2Cu308-5, the x and 
y axes become inequivalent, so that the s and d^2_j,2 
states on the one hand and the dxy x dx2-y2 and dxy 
states on the other are mixed. For La2-x Srx Cu O4 and 
Bi2Sr2CaCu208-i , it is the two orthogonal 45° axes which 
become inequivalent and then the mixing pairs are s and 
dxy onJjhe one hand and dx2-y2 and dxy x dx2-y2 on the 
other .eI 

Orthorhombic distortions on the square lattice should 
affect both the hopping integrals and the nearest neigh- 
bor interactions. These effects can be taken into ac- 
count by letting the hopping amplitudes ti and the 
nearest neighbor interaction Vi be different along the 
X and Y axes, as is the case of YBa2Cu308-5, or the 
45° axes (this would correspond to La2-xSrxCu04 and 
Bi2Sr2CaCu208-5). For example, for YBa2Cu308-5, the 
hopping integrals ti, t^, change to ij = ti[l ± a) 
a.ndVi^^x,y)=Vi{l±P). 

The expression for the current density 



2e 
ha 



dk 
(2^ 



^,<k + q\^k\' 



(17) 



has two terms. The main contribution comes from the 
condensate; the other is a back-flow term which comes 
from the quasiparticles and partially counteracts the con- 
tribution from the condensate. The quasi-particle term is 
also responsible for the non linear behavior at large values 
of q. The lattice spacing a has been set equal to the ex- 
perimental in-plane lattice constant of Bi2Sr2CaCu208-.5, 
which is 5.4 A. 
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The above expression for j* is valid to all orders in q 
and therefore goes well beyond linear response theory. To 
find the superfluid density tensor, j" should be expanded 
to linear order in q 

= -2ens,ij Vj,s, (18) 
so that 

ns,,=rn f l^d^ tanh (^) [M-% . (19) 

Here [Af^^] denotes the mass tensor. The current is 
therefore not parallel to the superfluid velocity, although 
the angle formed by them is always small. 

A mapping can be done from the tight-binding Hamil- 
tonian onto a much simpler linearized model £3 by ex- 
panding all functions and keeping only terms linear or 
quadratic in k. The following substitutions must then be 
made 

(-^ 2tjikp\ [k — kp)\ = fivp [k — kp) i— > £, 



400 




V„(meV) 

FIG. 1. Phase diagram for (a) the tight-binding model with 
doping x=0.2; (b) the linearized model, with parameters cho- 
sen to fit the order parameter fields along the axes. 

DISCUSSION OF ZERO CURRENT PROPERTIES 



r/a 



dk 



27r Jo 



^ A^{9) = v^{0) + y^{d) = Af + A^ cos(20) 



where A is some length scale of the order of the lattice 
constant. 

The saddle point equation for the order parameters 
fields can be written now as 



Af 



Jo 2^ Jo 



X tanh 



(3 {hvp ■ q+ Ve2 + |Ai(0)|2) \ 



(21 



where 
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24 



There are three parameters, namely kp X, g, and Ed, 
which can be adjusted to fit the values of A^ and A^ 
along the Vq and Vi axes. 



The numerical solution of the saddle points equations 
5b| ) for the order parameters is pretty straightforward. 
The zero temperature phase diagram is a function of 
the two coupling constants Vq and Vi, and the occupa- 
tion number n or, equivalently, the doping concentration 
(20) x = 1 — n. Fig. 1 (a) shows a representative case, taken at 
X = 0.2. There exist SDW, pure superconducting s aiid 
d, and mixed s + s* , s + s* + d, and d + i{s + s*) states.c^l 
The d-wave solution exists and is more stable close to 
half- filling, for values of Vi larger than | Vq I ■ The s-state 
solutions, on the contrary are more stable for small occu- 
pations of the band or when Vi < Vb ■ The time-reversal- 
symmetry broken solution d+i{s+s*) exists and is always 
more stable along the phase boundary between the d and 
the mixed s -I- s* states. The mixed s + s* + d only exists 
for large values of Vi and therefore doesn't appear in Fig. 
1(a). 

The effective on-site interaction in the case of the 
cuprates must still be repulsive and large even after a 
hypothetical process of renormalization of high frequency 
degrees of freedom. The part of the phase diagram rel- 
evant to these compounds is therefore the negative Vq 
quadrant, where the only physical phases are d-wave and 
SDW states. This model then nicely produces compe- 
tition between d-wave superconductivity and antiferro- 
magnetism, even though their phase boundary occurs for 
ratios Vi/|Vb| which are too large as compared with the 
"Mexican hat" effective potential referred to in the pre- 
vious section. The model also excludes any exotic time- 
reversal-symmetry broken state, like the d + i{s + s*). 
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FIG. 2. (a) Gaps Ad, s,s* ; (b) ratio 2Ati//('src as a function 
of doping X. 



Quantum fluctuations will likely melt the ordered phases 
into a strongly correlated electron liquid in the neigh- 
borhood of their phase boundary as in other ccftrelated 
low-dimensional systems with competing states.E3'ca In- 
deed, they could also displace the phase boundary to- 
wards larger values of |Vo|- The numerical results pre- 
sented for the tight-binding model in the rest of the paper 
will use the values (Vq, Vi) = (—150, 250) meV for the in- 
teraction constants, which correspond to a superconduct- 
ing solution in the neighborhood of the phase boundary 
between d-wave and SDW states. 

The doping dependence of the zero temperature gaps 
is plotted in Fig. 2 (a). The amplitude of the d-wave gap 
is always larger close at half-filling and vanishes at about 
x=0.7. The maximum occurs at x=0.22 irrespective of 
the values of the coupling constants due to fact that the 
Fermi level crosses over a van Hove singularity. Inter- 
estingly enough, there are s-wave solutions of the sad- 
dle point equations for large doping concentrations even 
though the on site interaction is repulsive. The s-wave 
gaps peak at a; ~ 0.8 because the Fermi level crosses a 
second van Hove singularity. 

The ratio 2 Ad/ifs Tc (Fig. 2 (b)) has a constant value 
of 4.3 for doping concentrations ranging from x=0.05 to 
0.4 (the ratio for a parabolic band is universal and equal 
to 4.25.) The highest value of Tc therefore occurs at 
x=:0.22, which is defined at the optimal doping concen- 
tration of this model. The ratio increases steeply for 
larger x, reaching a value of 7 when x is 0.67. Hence, the 
doping dependence of the critical temperaterie-jjesembles 
the universal curve Tc{x) of the cuprates.E3©E2l More- 
over, all physical magnitudes are expressed in real units, 
so that one may compare with experiments easily. The 



FIG. 3. (a) Normal pn{(^) and superconducting Ps{^) den- 
sities of states for three doping concentrations in a frequency 
range which covers almost all the band; (b) superconducting 
density of states ps(oj) for various x in the narrow frequency 
range [-200,200] meV. 



agreement doesn't extend too far though, as the exper- 
imental doping dependence of the zero tempep.ture gap 
is a monotonically decreasing function of x.O'Eil The ex- 
perimental ratio 2Arf(r — 0)/Kb Tc is therefore non con- 
stant. 

For underdoped cuprates there appears another char- 
acteristic temperature T*, which is a linearly decreas- 
ing function of the doping concentration. Above Tc but 
below T*, the d-wave gap partially closes down in such 
a way that the Fermi contour exists only in segments 
around the diagonals of the two-dimensional Brillouin 
zonclill Assuming that both T* and Ad{T — 0) vanish 
at a doping concentration of about 0.3, the dependence 
of both magnitudes with x is close to a straight line of 
slope pt,a ■ A possible generalization of 2 A /Kb Tc for the 
cuprates could therefore be the ratio 2*pa/pt, which is 
approximately constant and equal to about 5.5. 

The main source of trouble probably comes from the 
simple form of the dispersion relation of the model, where 
several drastic assumptions about the existence and na- 
ture of the hypothetical quasiparticles have been made: 
(a) There exists a well defined Fermi surface, even though 
this system has a pseudo gap above Tc. (b) Quasiparti- 
cles not only exist and are infinitely long-lived, but also 
their weight on the renormalized electronic propagator is 
always exactly 1. (c) The band shifts rigidly when the 
system is doped. To obtain the correct behavior for the 
doping dependence of Ad{T = 0), the treatment of the 
dispersion relation of the fermionic degrees of freedom 
should be improved. In particular, if one still wishes 
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-250 meV, 




25 


Vi = 150 meV and x 


= 0.2. 
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to use a quasi-particle picture, at least assumption (b) 
should be relaxed. 

The specific dispersion relation used gives rise to a 
second, weak van Hove singularity, in addition to that 
proper of the Hubbard t-t' model, see Fig. 3 (a). They are 
responsible for the maximums attained by the zero tem- 
perature gaps as a function of doping. The dispersion re- 
lation leads to a superconducting density of states which 
is largely asymmetric as is also seen in a-b tunneling-ep^ 
periments of Y Ba2 Cu3 Os-a and Bi2Sr2CaCu208-5J23'E2l 
Fig. 3 (b) is a plot of the density of states for three dop- 
ing concentrations close to the optimal value x=0.22, in 
the restricted frequency range which is often used in those 
tunneling experiments. 

The main effect of an orthorhombic distortion of the 
square lattice is to allow mixing of the d and s-wave 
states, as discussed above. Such an effect can indeed be 
seen in the numerical solution of the saddle point equa- 
tions, as shown in Table H. Another effects are related to 
the different role performed by the parameters a and f3. 
A finite a also reduces the value of the d-wave gap. On 
the other hand, a finite /3 doesn't alter the d-wave gap so 
much, but induces more appreciable s-wave components. 
The relative phase between s and d-wave gaps on two 
different twin domains can also be modeled by changing 
thgign of a and /3 as shown in the last four rows of table 

The mapping of the tight-binding model onto a lin- 
earized Hamiltonian works quite well. There are three 
parameters which appear naturally in the analytic calcu- 
lation, kp A, g, and Ed- All of them are needed to fit A^^ 
onto their partners, As^d along the axes. The values ex- 
tracted for kp A/tt, g, and Ed range from 1.21 toJJ27, 1.1 
to 1.6 eV-^, and 205 to 480 meV, respectively.^ There 
are no further parameters at hand so the mapping for the 
remaining part of the phase diagrams can be regarded as 
parameter free. Fig. 1 shows that the qualitative trends 
are indeed captured even though the slopes of the differ- 
ent phase boundaries are not exactly the same. 



INCLUSION OF SUPERFLOW 

The Ginsburg-Landau free energy functional of a pure 
d-wave superconductor on a perfect square lattice with a 
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FIG. 4. Zero Temperature gaps As,s*,d, free energy per site 
fsp, and current density j as a function of the wave vector q 
normalized to the size of the Brillouin zone for three different 
angles and doping concentrations (a) x=0.05; (b) x=0.2; (c) 
x=0.4. 
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is completely isotropic close to the critical temperature. 

1/2 

Here, the order parameter c>c [— |q!|/2/3] ' x 77^ x d 



1 /2 

a.nd S^s — [fi /2m I a I] is the superconducting coherence 
length. 

The relation between j and q is linear at the begin- 
ning, but for large enough q, j begins to saturate and 
attains a maximum which defines the critical current jc. 
Afterwards, the current density decreases and eventually 
vanishes, although this part of the curve corresponds to 
metastable states and has therefore no physical relevance. 
The linear behavior of j for small Doppler shifts arises 
from the dependence of the superfluid velocity with q, 
because Us is then almost constant. For Doppler shifts 
of the order of the inverse coherence length, on the con- 
trary, Hg is largely depressed and so the current density 
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FIG. 5. (a) Critical current density jc as a function of the 
angle 9q for three doping concentrations; (b) Angular spectral 
function A{u) — 0, for various values of q and 9q. Top panel, 
x=0.05, middle panel, x=0.2; bottom panel, x=0.4. 

saturates at a maximum which defines the critical current 
density jc- 

This qualitative behavior displayed by superconduc- 
tors close to Tc is also manifested at low temperatures. 
The top panels of Fig. 4 show the zero temperature gaps 
as a function of qa/2TT for three different angles 9q formed 
between g and the X axis, 0°, 22.5° and 45°. The most 
distinct feature is that the d-wave gap follows very sim- 
ilar trends for all the angles, showing that the superflow 
state is indeed rather isotropic. Besides, the gaps (a) 
decrease slowly as they would in the Ginsburg-Landau 
regime, but then drop abruptly at a large q value and, fi- 
nally, there is sometimes a strange tail; (b) small s and s* 
components are induced and grow slowly with q; (c) the 
range of g-values over which the gaps are finite is clearly 
larger for the nearly optimal doping concentration x=0.2; 
(d) the d-wave gap is somewhat smaller for 9q — 0° than 
for 9q — 22, 5°, 45°, but extends to larger values of q for 
x=0.05; on the contrary, it is bigger for the other two 
doping concentrations, but drops to zero earlier. 

The free energy is shown in the center panels. It is 
much the same for all angles at each doping concentra- 
tion and does not show any strange feature. The bot- 
tom panels show the current densities. They display the 
Ginsburg-Landau behavior in a loose sense, but there are 
again tails at the far right range of q. The overall form 
of the curves is fairly similar again for all the angles dis- 
played. Moreover, the superfluid density, which is the 
initial slope of the curve, seems to be the same for all 
of them. The form of the curves is markedly different 
for the different doping concentrations, and the current 



FIG. 6. Zero Temperature gaps As,s»,d, and current den- 
sity j as a function of the wave vector q normalized to the 
size of the Brillouin zone, for three possible values of the or- 
thorhombic parameters a and f5 and a doping concentration 
of x=0.2. The angles are (a) 6*, = 0° and (b) = 45°. 

densities attained are larger again for x=0.2. 

The critical currents depend only slightly on the an- 
gle 0q, as shown in Fig. 5 (a) for several doping concen- 
trations. A closer look reveals that jc is maximum at 
dq = 45° and minimum at 9q — 0° for x=0.05; it is still 
maximum at 45 degrees at x=0.2, but in this case jc is 
a non-monotonic function of the angle; the maximum is 
attained at 9q = 0° for x=0.4. This behavior can be 
understood in terms of the angular spectral function 

A(c^,C<45°)= / k5{Lo~^{kX)), (24) 

which weights the back-flow term due to quasiparticles 
in Eq. (p^). A{uj — 0, C) (Fig. 5 (b)) has strong peaks for 
small angles when x=0.05 and 0.2 and therefore the back- 
flow contribution to the current must then be important. 
A{u} = 0, C) is smoother for x=0.4, and accordingly, jc 
is flatter. Indeed, the angular dependence of the critical 
current jc is even weaker for the linearized model; jc is 
always minimum at 9q — 45° in such a case. Taken to- 
gether, all these facts clearly indicate that the anisotropy 
of the band is a more efficient source of anisotropy than 
the internal structure of the order parameter. 

The critical currents achieved are of the order of 10 
Amp/cm. This leads to an estimated critical current 
for Bi2Sr2CaCu208-5 of 5 x 10^ Amp/cm^, which is one 
order of magnitude larger than tke. value 10^ Amp/cm^ 
found in thin films some time ago.EZI Similarly, we obtain 
jc — 2 X 10® Amp/cm^ for YBa2Cu308-5, while the ex- 
perimental values are of the order of 5 x 10^ Amp/cm^, 
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not far away.EzI 

The superfluid density tensor 



^|„,|na.l,(iS)[M-.]„ ,25) 



is not proportional to |Ap as happens in Ginsburg- 
Landau theory. 

Orthorhombic effects modeled by the /3 parameter do 
not affect the behavior of the d-wave gap or the current 
densities. They only induce a rather large s component, 
which also increases with q. The curves for A^; and j are 
indeed indistinguishable from those at /3 = 0. The impact 
of the other orthorhombic parameter a is bigger. Fig. 6 
shows the gaps and current densities for two opposite or- 
thorhombic elongations of the X and Y axes, a = ±0.05. 
For an angle 6q = 0°, (a) the d-wave gaps are initially 
reduced by the same magnitude as they should, but as q 
increases the solution for negative a drops to zero faster; 
(b) the s-wave gaps have initially the same magnitude 
and opposite signs, but that corresponding to a = —0.05 
goes through zero, changes its sign and eventually drops 
to zero again; (c) the current densities show the oppo- 
site trend, in the sense that both the superfluid density 
Us and the critical current jc are larger for negative a. 
Moreover, superfluid densities and critical currents for 
a = 0.05 are even smaller than those corresponding to 
no orthorhombic distortion. 

Fig. 6 (b) corresponds to 9q = 45*^. In this case, the two 
orthorhombic distortions corresponding to a = ±0.05 are 
equivalent and the curves for the d-wave gaps fall one on 
top of the other. The two s-wave gaps have exactly the 
same magnitude and different sign. The current densities 
are also equal as they should. The overall trend is that 
a finite a reduces both the d-wave gap and the critical 
current, although the superfluid density does not change 
much. 

Comparison of Figs. 6 (a) and (b) shows that or- 
thorhombicity is also a more efficient source of anisotropy 
that the internal symmetry of the order parameter func- 
tion. Indeed, not only jc but also ris strongly depend on 
9q even for these small values of the parameter a. 

The time- reversal-symmetry broken solution d + i{s + 
s*) to the saddle point equations appears for attractive 
Vo around the phase boundary between the d and s + s* 
solutions. In this part of the phase diagram, there exist 
not only the d+i{s + s*) state, but also flnite s + s* and d 
solutions. The phase transitions from one stable state to 
another are therefore first order. A finite Doppler shift 
gives rise to (a) mixing between the d and s + s* states so 
that each of these two solutions is of the form s+s* + [d] or 
[s+s*]+d; where the subdominant component is enclosed 
in brackets; (b) competition between these two solutions 
and with the broken time-reversal-symmetry state. 

Numerical calculations for the linearized case and suit- 
ably chosen values of Vj^ and Vj" are displayed in Figs. 
7 and 8 for angles 6q = O" and 45°. They show that 



0.008 

0.004 

0.000 

-0.004 

-0.008 
0.00 



-0.02 



^ -0.04 



-0.06 

•0.2 
s 

0.1 

■s 

0.0 



(a) 



s part of [d]+s 
d part of [d]+s 
s part of d+[s 
d part of d+is 
s part of d+[s] 
d part of d+[s] 



[d]+s 
d+is 
d+[s] 



(C) 



[d]+s 
d+is 
d+[s] 



0.000 0.004 0.008 

^sP/Ed 



0.012 



FIG. 7. (a) Zero temperature gaps Ai, (b) free energy fsp, 
and (c) current density j in arbitrary units as a function of 
the normalized superfluid velocity in the linearized model for 
= 0° and g^ = 0.15, Vj- g^ = 0.32; d + [s] and [d] + s 
indicate mixed states with dominant d or s components, re- 
spectively. 



the d + is state is initially more stable but the s com- 
ponent decreases when Vs increases and eventually van- 
ishes. The d + is state disappears at that moment in 
favor of the d + [s] state, [d] + s being higher in energy. 
The d + [s] solution becomes metastable once the critical 
current density of that state is reached. There is then a 
phase transition to the [d] + s state, which has a higher jc, 
because it is the current what is actually imposed in an 
experimental setup. Above this second critical current 
density, the system reverts to its normal phase. 



CONCLUSIONS 

The formalism needed to include superflow in a gauge 
invariant way in a superconductor has been discussed 
in detail and shown to be equivalent to give a Doppler 
shift to the wave vectors in the kinetic energy part of the 
Hamiltonian. This has been applied to study supercur- 
rents in a clean, infinite, d-wave superconductor on an 
orthorhombically distorted square lattice, using a model 
which describes some of the features of the superconduct- 
ing state of several High-temperature superconductors. 
The behavior of gaps, free energy, current density and 
critical current as a function of modulus and angle of the 
Doppler shift, doping, and orthorhombic distortion have 
been discussed. 

There are three possible sources of anisotropy in the 
discussed model: the internal symmetry of the order pa- 



9 



0.006 



s part of [d]+s 
d part of [d]+s 
s part of d+is 
d part of d+is 
s part of d+[s] 
d part of d+[s] 




0.000 



0.002 0.004 



0.006 



0.008 



FIG. 8. (a) Zero temperature gaps A;, (b) free energy fsp 
,and (c) current density j in arbitrary units as a function of 
the normalized superfluid velocity in the linearized model for 
6»q = 45° and Vj^ = 0.15, Va g^ = 0.32. 



rameter, the complex dispersion relation of the band and 
the orthorhombic distortion of the lattice. The results 
of this article show that the most important one is or- 
thorhombicity distortion, then the form of the band and 
last and least the d-wave form of the gap function. 

After this work was completed wc foimd that D. Fcdcr 
and C. Kallin had also been studying the behavior of su- 
percurrents in d-wave superconductors and had reached 
conclusions which agree with those of this article. 
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